% Copyright (C) 2021-2023 Benjamin Born, Francesco D'Ascanio, Gernot J. Mueller, Johannes Pfeifer
%
% This is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% It is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
% 
% For a copy of the GNU General Public License,
% see <http://www.gnu.org/licenses/>.


transp90 = 0.25;
linW = 2;

FX_pos=4;
G_pos=2;
%flip sign to be consistent with notation
IRF_dist.BP_simple=quantile(IRF.BP_simple,ptileVEC,4);
IRF_dist.BP_simple(:,FX_pos,:,:)=-IRF_dist.BP_simple(:,FX_pos,:,:);
IRF_dist.PF_simple=quantile(IRF.PF_simple,ptileVEC,4);
IRF_dist.PF_simple(:,FX_pos,:,:)=-IRF_dist.PF_simple(:,FX_pos,:,:);

var_names={'T','G','Y','FX (down is apprecation)'};

method_names={'BP_simple','PF_simple'};
shock_names={'T shock','G shock'};
for method_iter=1:length(method_names)
    hh=figure('Name',method_names{method_iter});
    for shock_iter=1:2
        for var_iter=1:n
            subplot(6,2,var_iter*2-(shock_iter-1))
            hline(0,'k');
            hold on
            a = squeeze(IRF_dist.(method_names{method_iter})(:,var_iter,shock_iter,1))';
            b = squeeze(IRF_dist.(method_names{method_iter})(:,var_iter,shock_iter,end))';
            [~,~]=jbfill(0:1:Horizon-1,a,b,'b','none',1,transp90);
            hold on
            plot(0:1:Horizon-1,squeeze(IRF_dist.(method_names{method_iter})(:,var_iter,shock_iter,3)),'b','LineWidth',linW);

            if shock_iter==2 && (strcmp(method_names{method_iter},'BP_simple') || strcmp(method_names{method_iter},'PF_simple'))
                plot(0:1:Horizon-1,squeeze(IRF_dist.(method_names{method_iter})(:,var_iter,3,3)),'r','LineWidth',linW);
            end
            set(gca,'XTick',0:4:Horizon-1)
            xlim([0 40])
            box on
            if var_iter==1
                title([shock_names{shock_iter} ': ' var_names{var_iter}])
            else
                title(var_names{var_iter})
            end
        end
    end
    dim = [12,8];
    set(hh,'paperpositionmode','manual','paperunits','inches');
    set(hh,'papersize',dim,'paperposition',[0,0,dim]);
    print(hh,'-dpdf' ,['./Figures/' method_names{method_iter}]);
end

G_shock_pos=2;
scaling_factor_BP=IRF_dist.BP_simple(1,G_pos,G_shock_pos,2);
scaling_factor_PF=IRF_dist.PF_simple(1,G_pos,G_shock_pos,2);
BP_G_IRFS=squeeze(IRF_dist.BP_simple(:,[2,3,4,1],G_shock_pos,:))/scaling_factor_BP;
PF_G_IRFS=squeeze(IRF_dist.PF_simple(:,[2,3,4,1],G_shock_pos,:))/scaling_factor_PF;
var_names=var_names([2,3,4,1]);